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Abstract 

Junctions  are  important  features  for  image  analysis  and  form  a  critical  aspect  of  image  under¬ 
standing  tasks  such  as  object  recognition.  We  present  a  unified  approach  to  detecting  (location  of 
the  center  of  the  junction),  classifying  (by  the  number  of  wedges  -  lines,  corners,  3-junctions  such 
as  T  or  Y  junctions,  or  4-junctions  such  as  A-junctions)  and  reconstructing  junctions  (in  terms  of 
radius  size,  the  angles  of  each  wedge  and  the  intensity  in  each  of  the  wedges)  in  images.  Our  main 
contribution  is  a  modeling  of  the  junction  which  is  complex  enough  to  handle  all  these  issues  and 
yet  simple  enough  to  admit  an  effective  dynamic  programming  solution.  Broadly,  we  use  a  template 
deformation  framework  along  with  a  gradient  criterium  to  detect  radial  partitions  of  the  template. 
We  use  the  minimum  description  length  principle  to  obtain  the  optimal  number  of  partitions  that 
best  describes  the  junction.  Kona  [27]  is  an  implementation  of  this  model.  We  (quantitatively) 
demonstrate  the  stability  and  robustness  of  the  detector  by  analyzing  its  behavior  in  the  presence 
of  noise,  using  synthetic/controlled  apparatus.  We  also  present  a  qualitative  study  of  its  behavior 
on  real  images. 


I.  Introduction 

A  critical  component  of  most  recognition  systems  is  a  stable,  representative  feature  extraction  from 
images.  One  of  the  key  features  used  in  recognition  is  junctions:  T-junctions,  Y-j unctions,  Ad¬ 
junctions  and  so  on.  These  junctions  are  also  critical  for  stereo  vision  modules  or  motion  modules, 
since  these  are  places  where  occlusions  can  be  identified.  Such  points,  for  example,  coincide  with  the 
images  of  trihedral  vertices  of  an  object.  These  are  critical  features  for  recognition  as  suggested  by 
[6],  [7],  [33].  They  are  also  critical  features  for  motion  as  suggested  by  Walach  experiments  [34]  and 
by  Movshon  and  Adelson  [1].  A  study  of  the  role  of  junctions  in  stereo  is  presented  in  Malik  [21]. 

There  have  been  basically  two  different  views/approaches  for  detecting  junctions:  edge  detection 
followed  by  grouping  of  edges  to  form  junctions  [25],  [24],  [4],  [2],  [13],  and,  treating  a  junction 
as  a  template  matching  phenomenon  [8],  [12],  [15].  In  the  former,  it  is  assumed  that  the  presence 
(or  absence)  of  a  junction  is  determined  by  “grouping”  the  intensity  gradients  near  a  hypothesized 
junction.  Usually  one  is  interested  in  examining  large  gradients  in  the  direction  perpendicular  to 
the  hypothesized  Radial  line.  Experiments  in  this  framework  are  limited  and  even  the  richest  ones 
shown  in  [25]  are  interesting  but  not  exhaustive.  In  the  latter  approach  it  is  assumed  that  a  (suitably 
small)  local  neighborhood  is  sufficient  to  detect  a  junction*  The  basic  idea  is  to  fit  a  junction-model 
to  the  input  signal  in  a  neighborhood.  This  involves  minimizing  an  energy  function  which  gives  a 
measure  of  the  “distance”  of  the  junction-model  from  the  input  signal. 
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We  use  the  template  deformation  framework  to  develop  a  “junction  detector,”  to  find  corners  (two- 
junctions),  tri-corners  (tri-junctions),  quad-corners  (quad-junctions),  etc.,  defined  as  points  where 
two  or  more  homogeneous  surface  patches  are  located  within  a  neighborhood  of  the  point.  Our  ap¬ 
proach  is  to  use  a  combination  of  the  two  paradigms:  grouping  of  edges  (via  Dynamic  Programming) 
and  fitting  templates.  We  use  a  template  deformation  framework,  using  the  minimum  description 
length  (MDL)  principle,  and  we  include  the  gradient  criterium  in  order  to  detect  the  radial  partitions 
of  the  template  as  a  grouping  mechanism.  The  task  is  to  find  the  minimum  number  of  wedges  that 
best  describes  the  junction.  Note  that  as  we  increase  the  number  of  wedges,  the  junction  description 
gets  more  accurate;  hence  the  task  is  to  use  the  MDL  principle  to  obtain  an  optimal  number  of  such 
wedges.  The  minimum  length  encoding  principle  was  first  suggested  by  J.  Rissanen  [28],  [30],  [31], 
and  is  very  close  in  spirit  to  the  Kolmogorov  sufficient  statistic,  as  discussed  in  [3]. 

Kona  1  [27]  is  an  implementation  of  our  model.  To  test  the  stability/robustness  of  the  detector  we 
have  experimentally  analyzed  its  behavior  for  the  location  and  junction  parameters  against  noise.  We 
used  a  synthetic/controlled  apparatus  for  a  quantitative  study.  We  have  also  qualitatively  studied 
the  stability/robustness  of  detector  for  real  stereo  images. 

Our  model  is  in  contrast  to  other  approaches  which  we  briefly  review  here.  In  [13]  ,  in  order  to  find 
corners,  a  local  operator  is  first  applied  to  detect  possible  junction  locations.  Then,  inside  an  area 
around  this,  the  edge  lines  are  detected  and  depending  on  this  number  (of  edge  lines),  a  parametric 
model  is  used  to  fit  the  data  is  fit  in  a  small  region  (25  X  25  pixels). 

The  idea  of  performing  local  feature  detection  by  projecting  image  data  onto  a  subspace  is  funda¬ 
mental  in  [8],  [10].  Basically,  the  input  is  orthogonally  projected  onto  a  finite  dimensional  subspace 
of  the  Hilbert  space  of  functions.  An  energy  function  (which  is  the  L 2  norm  of  the  the  difference  of 
the  input  and  the  fitted  function)  is  minimized  in  this  finite  dimensional  space.  The  two  main  issues 
are  finding  an  orthonormal  basis  that  spans  a  good  finite  dimensional  subspace  and  minimizing  the 
energy  function.  This  approach  can  give  closed  form  solutions  for  edges  [8],  [10],  and  lines  [9].  A 
generalization  to  junctions  is  suggested  in  [32],  using  steerable  filters,  and  in  [18],  through  the  use  of 
principal  component  analysis  (PCA).  A  possible  drawback  for  the  PCA  is  that  it  does  not  have  an 
explicit  model  of  the  junctions  (to  be  adapted). 

The  authors  in  [15],  [12]  employ  an  iconic  model  of  the  junctions,  giving  forms  similar  to  our  1Z 
(equation  (5)  in  Section  II).  See  Section  II  for  a  comparison. 

In  [16],  corners  and  junctions  (which  are  modeled  as  two  adjacent  corners)  are  represented  by 
functions  (models)  that  are  blurred  with  a  Gaussian  (or  an  exponential  filter)  where  the  authors  use 
a  closed  form  solution.  In  general,  numerical  methods  are  used  to  obtain  parameters  that  minimize 
the  distance  to  the  input  data  using  an  L 2  norm.  This  is  also  the  case  in  [17]. 

The  paper  is  organized  as  follows:  Section  II  describes  the  junction  model.  Section  III  studies 
detailed  issues  of  the  optimization  process  to  estimate  the  junction  parameters.  Section  IV  presents 
the  results  with  a  study  of  the  stability  of  the  detector.  Section  V  concludes  the  paper. 


^he  word  for  “corner”  in  Hindi  is  Kona. 
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II.  The  Junction  Model 

We  model  a  junction  as  a  region  of  an  image  where  the  values  are  piecewise  constant  in  wedge- 
shaped  regions  emanating  radially  from  a  central  point,  covering  a  small  disk  centered  at  the  point 
and  omitting  a  (much)  smaller  disc  centered  at  this  point  (see  Figure  1).  The  parameters  of  a  junction 
consist  of  (i)  the  radius  of  the  junction-disk,  (ii)  the  center  location,  (iii)  the  number  of  radial  line 
boundaries,  (iv)  the  angular  direction  of  each  such  boundary,  and  (v)  the  intensity  within  each  wedge. 
The  radius  of  the  disk  addresses  the  “scale”  issue,  and  the  location  of  the  center  is  a  kind  of  “interest 
operator”  [11]  that  determines  the  position  where  the  feature  is  located  in  a  region,  possibly  pre¬ 
defined.  Another  corner-detector  by  [35]  performs  a  similar  role  of  detecting  locations  that  have  a 
corner.  This  uses  a  predetermined  window/mask  on  each  pixel  and  compares  the  intensity  of  every 
pixel  with  that  of  the  center  of  the  mask  to  determine  whether  there  is  a  corner  at  the  location. 

We  can  formulate  the  junction  detection  problem  as  one  of  finding  the  parameter  values  that  yield 
a  junction  that  best  approximates  the  local  data  using  minimum  description,  and  declaring  local 
minima  of  the  error  as  junctions.  The  best-fit  parameter  values  provide  attributes  of  the  detected 
junction. 

Let  T  denote  the  piecewise  constant  function/template.  It  has  N  angles  and  N  intensities  if  N  is 
the  number  of  constant  pieces.  Further,  let  /  denote  the  input  signal. 

Define  the  energy  function,  at  a  point  (i,  j)  on  the  image  as  follows 

E  =  V  +  \Q ,  (1) 


where  A  >  0. 

The  first  term,  V1  is  a  measure  of  the  distance  of  the  fitted  function  from  the  data  using  the  L 2  norm: 

PC O  P  27T 

V=  /  \I(r,8)-T(6)\2g(r)rdrd0,  (2) 

Jo  Jo 

where  g(r)  is  an  appropriate  modulating  function  that  goes  to  zero  for  large  r,  thus  defining  the 
template  size. 

We  note  that  the  input  to  the  system  is  an  image  Iraw  that  can  be  filtered  by  a  Gaussian  mask  to 
yield  a  smooth  image  I.  Smoothing  Iraw  is  needed  due  to  the  discretization  of  the  image  lattice  and 
how  it  affects  the  junction  detection.  Our  analysis  of  this  effect  is  shown  in  Section  IV-B. 

The  second  term,  Q)  is  a  measure  of  the  distance  of  the  gradient  using  the  L2  norm. 

roo  p2it 

Q=  /  /  |V/(r,0)-  VT(8)\2g*(r)rdrd8,  (3) 

Jo  Jo 

where  g*(r)  is  an  appropriate  modulating  function,  not  necessarily  the  same  as  g(r). 
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v/  = 


dl(r,8)  ,  19/M) 

~^~Gr  +  r  1)8  ~e° 


Note  that 
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where  er  and  e#  are  the  orthonormal  vectors  in  the  r  and  9  direction  respectively,  evaluated  at  (r,  0). 
Separating  the  angular  and  the  radial  terms,  Q  can  be  written  as, 

Q  =  Aed  +  7 £er, 

where, 


(4) 


/•co  f  2tt  f  QI\2 

1Z  =  L  l  (w)  (5) 

Taking  </*(r)  =  r2G<j(r)  where  G( T(r)  is  the  Gaussian  function  with  standard  deviation  <j,  our  function 
7v()  is  the  same  as  the  regularity  measure  5()  for  junctions  as  in  equation  (4)  of  [15].  Also,  it  is  the 
same  as  UJ&  defined  in  equation  (31)  in  [12]  which  is  the  measure  of  spirality  omitting  the  term  normal 
to  the  radial  direction.  However,  in  [12],  a  general  method  is  presented  to  detect  circular  symmetry; 
thus  our  function  1Z  could  be  considered  as  a  special  case. 


A.  On  the  Exact  form  of  Template  T  &  Energy  E 

This  section  gives  a  more  rigorous  definition  of  T  and  E  as  a  continuous  function  of  two  variables. 
The  reader  may  skip  this  section  without  loss  of  continuity. 

Let  the  image  be  given  by  I{x,y).  T  is  a  template  that  partitions  the  image  in  N  regions,  an  “In¬ 
junction”  template.  Between  the  partition  lines  (wedges),  say  p  and  p+ 1,  we  assign  a  constant  value 
Tp  for  the  template.  Let  us  further  consider  a  template  for  the  edges  that  is  zero  everywhere  except 
along  the  partition  lines,  6  =  6p.  We  can  write  both  templates  as 

N 

f(»)  =  {{Tv,6vy)P£(l,2,...,N)}  =  YJTvu(6-0v)u{6v+1-0) 

P=  1 
N 

T'(0)  =  e(e)  =  {(ep,0py,pe  (1,2, . . N)}  =  *£  ep5(0  -  6P) ,  (6) 

P=  1 


where  0jv+i  —  #i,  S(x)  is  the  Dirac  delta  function  and  u(x)  is  the  step  function  (or  Heaviside 
function),  i.e. 


/oo 

<5(:r)  dx  =  1  <5(ic) 

-oo 


co  x  =  0 

0  X 


fx 

and  u(x)  —  /  S(x)  dx  — 

J  —  CO 


1  X  >  0 

0  x  <  0 


To  fit  a  junction  template  T(B)  to  the  image  I(x ,  y)  and  to  fit  an  edge  template  e(6)  to  the  image 
gradient  along  the  angle  6  (perpendicular  to  the  radial  line),  an  error  function  based  on  the  L2  norm 
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weighted  by  a  function  g(r)  is  used.  Using  the  polar  coordinate  system, 

r  =  J  x2  +  y 2  ,  6  =  arctan  —  ,  and  x  =  r  cos  6  ,  y  =  r  sin  6  , 
v  x 

the  image  gradient  is 

0J(r,0)  1  dl(r,0) 

Vi  =  — - er  H - — — eg. 


dr 


r  80 


Thus  the  energy  E  at  a  point  (i,j}  is: 


E  =  6)- T{e)fg(r)  +  \\E-E^-e(e)\^g*{r)}rdrde,  (7) 

where  g(r)  and  g*(r)  are  functions  that  weights  the  importance  of  the  piecewise  constant  fit  according 
to  the  distance  from  the  center  r  =  0. 


III.  Energy  Minimization 


Recall  the  energy  equation: 


E  =  A  U+(V  +  XA) 

=  XR  +  E.  (8) 

1Z  is  independent  of  the  junction  template  and  is  used  to  filter  out  the  image  locations.  E  is  minimized 
to  obtain  the  most  appropriate  junction  parameters.  We  now  discuss  the  details  of  computing  1Z  (scale 
and  location)  and  E  (junction  parameters). 


A.  Scale  and  Location  (OnlZ) 


We  have  considered  g(r)  =  1/r.  Notice  that  for  A  to  be  used  in  equation  (1),  g*(r)  =  r2g(r)  =  r  is  a 
choice  that  makes  V  and  Q  have  the  same  unit  (intensity2/length).  (Thus  A  is  unitless  as  expected.) 


Estimating  Ro  and  R\  (scale):  We  further  refine  g(r)  as  follows: 

g(r)  = 


0 

1 


0 


r  <  Ro 
Ro  <  r  <  R\ 

r  >  Ri 


Notice  that  we  introduce  a  “hole”  of  size  Rq  >  0  and  R\  is  the  size  of  the  window.  To  compute  Rq 
and  Ri,  we  study  7\U,  which  is  1Z  for  r,  the  radius  being  considered.  Figure  2  shows  the  plot  of  the 
relative  values  of  7^r,  er,  which  is 

7?r 

(9) 


e  = 


nr+1' 

Notice  that  7 Z  has  comparatively  larger  values  close  to  the  origin,  which  we  remove  by  introducing  a 
“hole”.  We  look  for  a  range  of  r  values  where  er  is  less  than  a  threshold,  rw.  This  leads  to  obtaining 
two  values  that  correspond  to  Rq  and  R\.  See  Figure  2  and  3  for  examples. 


The  weight  function  a(x)  which  is  defined  in  equation  (37)  in  [12]  is  similar  to  our  g(r)  function. 
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The  use  of  a  hole  is  also  supported  by  psychophysical  experiments  [34]  that  suggest  humans  may  also 
utilize  such  a  method.  The  experiments  suggested  by  Wallach  [34]  have  shown  that  different  motion 
effects  can  be  produced  depending  on  the  junction  movement.  In  a  Wallach  type  of  experiment, 
junctions  occur  due  to  the  intersection  of  lines,  and,  junctions  move  due  to  the  movement  of  the 
lines.  The  results  of  the  experiments  in  [34]  seem  to  suggest  that  junctions  are  detected  although 
there  is  a  small  gap  (Ro)  between  the  lines. 

Estimating  the  location:  We  do  the  following  to  select  the  best  location  in  an  image  region2.  We 
compute  W  =  IZ/r  for  all  the  points  in  the  region  (with  not  necessarily  the  same  window  size).  The 
one  with  the  minimum  value  defines  the  location.  We  illustrate  this  in  Figure  4. 


B.  Reconstructing  the  junction  (On  E) 


In  this  section  we  explain  how  the  junction  parameters  (number  of  wedges,  wedge  angles  and  wedge 
intensities)  are  estimated.  Recall  that  the  last  section  discussed  the  estimation  of  the  radius  (size) 
of  the  junction  and  the  size  of  the  hole. 


We  can  carry  out  an  appropriate  numerical  integration  to  obtain  the  value  of  the  first  term  7 Z.  The 
second  factor,  if,  is  used  to  estimate  the  junction  template.  The  unknown  factors  are:  TV,  the  number 
of  intersecting  lines  (or  wedges)  at  the  junction  3,  {0p},  {Tp},p  =  1 .  . .  TV,  where  N  is  the  number  of 
wedges,  6p  s  are  the  angles  where  the  partitions  occurs,  Tp  s  are  the  intensity  values. 


We  can  write  down  E  as 


E  =  T  +  V, 


where,  T  is  fixed,  that  is,  it  does  not  depend  on  the  unknown  parameters,  whereas  V  does.  Recall 
g*  (r )  =  r2g(r).  With  some  straightforward  manipulations  we  obtain 


T  = 


OO  /*27T 


0  JO 


/2(r,  0)  +  A( 


dl(r,8) 

d6 


V  = 


poo 

p  2tt 

L 

L 

T2(6)  -  2T  (6)1(6)  +  A 


8T 


81  8T 


86  2  86  '86 


g(r)rdrd6 

g(r)rdrd6 


T  can  be  approximated  numerically.  Also,  V  can  be  approximated  as  V. 

N 

v  =  £  [  (%,  -  ep)(-2TvhpMl  +  T2C)  +  A(-2 T‘v(dei)e„  +  C*+) 

P=  1 


where, 


poo 

C  =  /  g(r)rdr  =  Ri  -  R0,  (10) 

Jo 

poo 

1(0)  =  /  l(r,0)g(r)rdr,  (11) 

Jo 


2Notice  that  7Z  has  unit,  (intensity)2 length;  thus  7Zr  has  unit  (intensity)2. 

3N  for  homogeneous  region  is  1,  for  line  and  corner  is  2,  for  junctions  like  T-,  Y-junction  is  3,  4  for  X-junction  and 


so  on. 
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POO  POO 

C*  =  I  g*(r)rdr=  /  r2  d 
Jo  Jo 


r  = 


p3  p3 
/7l  U0 


mh  =  d-M 

{9e  )ep  -  dfj 


9JT  7(r.  %(r)rdr 

de 


h+ 1 

'5,  A-:.  =  A0£/(iA0).  (12) 

3  ~jp 

Further,  7(0),  (<9^I)gF  can  be  approximated  numerically.  For  the  sake  of  brevity,  we  have  omitted  the 
details  of  the  derivations  here,  and  present  it  in  Appendix  A. 


Estimating  wedge  angles  and  intensities:  Although  the  energy  equation  E  looks  fairly  complex,  it  has 
a  remarkably  simple  and  natural  interpretation.  Once  the  role  of  each  factor  of  the  energy  equation 
is  ascertained,  some  computations  (like  evaluating  T ,  for  instance),  can  be  omitted  without  changing 
the  solution.  To  recapitulate,  the  factors  are  used  thus: 

1.  1Z  is  used  for  scale,  i.e. ,  to  determine  the  size  of  the  window,  R\.  This  is  also  used  to  obtain 
the  exact  location  of  the  junction  in  a  neighborhood. 

2.  E  is  used  to  obtain  the  junction  parameters. 

The  1(6)  (equation  (11))  can  be  viewed  as  integrating  the  intensity  along  a  radial  line.  Thus,  the 
two-dimensional  image  is  projected  on  to  a  one-dimensional  coordinate  0,  appropriately  discretized. 

Let  I(6i)  be  defined  for  02l  . .  . ,  0^  and,  (6t+1  -  0J  =  2n fd,  Vi  4.  For  a  p-junction,  I\,  T2, . .  . ,  Tp, 
are  the  template  intensities  and  the  wedge  boundaries  are  at  0^  ,  0&2 *  •  •  • ,  @kpi  Ofa  <  0fc2  ■  ■  ■  <  @kp- 

Assuming  we  know  the  0p’s,  we  can  obtain  Tp s  by  setting  =  0,  Vp.  When  A  =  0,  —  f(6p ), 

for  some  function  /().  Thus  fixing  6p  s,  the  Tf' s  can  be  shown  to  be  the  following: 

T  g UfAf 

k+i  -  h 

In  other  words,  T\  is  a  piecewise  constant  fit  which  is  the  average  value  of  the  data  in  that  region. 
The  energy  for  the  fit  is 

P  kl  +  1  /  x2 

r = e  E  {m  -  /;)  • 

/=!  j^ki 


When  A  ^  0,  ~E?  has  some  extra  terms.  It  can  be  verified  that  T(  is  the  same  as  before,  and,  EF  is: 

^  =  tf  (h«4  -  tf + A  E  (h«4  -  7W-i))2- 

/=!  j—ki  i—k\ 


We  compute  the  6p  s  by  exploring  all  possible  set  of  6p  s.  We  summarize  the  dynamic  program 


cost  of  fitting  Tp  to  6j1  0j+ 1  ...  0/,  l  >  j, 

cost  of  fitting  Tp  to  0/,  0j+1  .  .  .  6\  . .  .  0y ,  l  <  j. 


4  If  *  =  d,  i  +  1  is  defined  to  be  1. 
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f 

P—1, 

\  mhl?</ 

V*<U>) 

otherwise. 

E  =  min  £™m 

s—l,m 

This  dynamic  programming  solution  has  been  implemented  in  Kona. 

T  is  a  constant  (fixed),  depending  only  on  the  image,  and  does  not  influence  the  junction  param¬ 
eters.  The  scenario  where  T  could  play  a  role  is  at  the  time  of  comparing  the  energy  values  in  a 
neighborhood.  But  we  use  1Z  to  filter  the  neighborhood  so  T  may  be  ignored  without  any  damage. 

Estimating  the  number  of  wedges:  We  estimate  the  number  of  wedges  by  measuring  the  rate  of  increase 
of  energy  by  thresholding  the  ratio  of  the  energy  measurement,  En .  More  precisely,  the  (optimal) 
number  of  wedges,  V,  is  computed  by  thresholding  the  relative  error,  rn, 

~Wn+1 

r"  =  C^. 

jp 

Notice  that  by  definition  rn  is  unitless.  Although,  in  principle,  we  are  looking  for  the  minimum  rn, 
in  practice  we  terminate  the  computation  when  rn  drops  below  a  pre-defined  threshold,  say  rc,  as  n 
is  increased.  (See  Figure  5  for  an  example.)  Note  that  as  the  number  of  parameters  increases,  ji'n 
decreases,  i.e. ,  rn  <  1. 

To  summarize,  the  following  steps  are  involved  in  detecting  junctions  on  a  large  region  of  an  image: 

1.  Compute  7£,  the  measure  of  radial  variation,  and  Ro,  the  size  of  the  hole  and  i?i,  the  size  of  the 
template  at  every  point. 

2.  Filter  the  locations  using  a  threshold  on  lZr . 

3.  In  a  neighborhood  of  a  filtered  location,  pick  the  one  with  minimum  7 7,  and  remove  all  other 
locations  within  a  radius  of  Ri  of  this.  Repeat  this  for  all  the  filtered  locations. 

4.  Compute  the  junction  parameter  for  all  the  filtered  locations. 

Kona  is  an  implementation  of  this  junction  model.  It  is  programmed  in  C.  A  version  in  C++  and  a 
parallel  implementation  of  Kona  is  also  available. 

IV.  Results  of  Experiments 

We  carried  out  a  series  of  experiments  on  synthetically  generated  images  (see  Figures  6  and  7)  to 
test  the  stability  of  the  algorithm  in  the  presence  of  noise  in  the  image,  and  to  further  understand 
the  role  of  smoothing.  In  these  experiments,  we  suppressed  the  use  of  threshold  values  using  a  fixed 
Ri  =  4  and  fixed  N  =  3,  i.e.,  forcing  the  algorithm  to  pick  up  the  optimal  3-corners.  We  then  carried 
the  experiments  on  real  images  to  test  the  general  performance  of  the  algorithm. 


A.  Stability  of  the  Algorithm 

Figure  6  shows  the  result  of  one  such  experiment  where  we  look  for  3-corners  in  the  center  region  of 
the  image.  It  shows  that  when  the  image  edges  are  sharp  (<j  =  0  .  .  .3  in  the  figure),  the  error  in  the 
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angles  and  the  intensities  are  slightly  higher.  The  errors  are  least  in  the  range  of  a  =  6  . .  .20,  and 
increase  further  on. 

Note  that  the  (Manhattan)  distance  of  the  location  of  the  3-corner  is  less  than  6  pixels;  the  error 
in  the  angles  is  bounded  by  30°  at  the  very  worst  and  the  intensity  differs  by  40  units  in  the  worst 
situation.  It  also  shows  that,  in  the  best  situation  the  angle  differs,  on  an  average,  by  15°  from  the 
true  answer. 

B.  Effect  of  Image  Smoothing 

We  present  in  Figure  7  the  results  of  one  such  experiment.  The  images  used  here  are  the  smoothed 
version  of  the  ones  shown  in  Figure  6.  Note  the  reduction  in  error  around  the  sharp  images  or  low 
values  of  u. 

It  is  worth  noting  that  even  for  sharp  junctions  (without  noise  added)  the  algorithm  performs  better 
when  smoothing  is  introduced.  It  is  possible  that  smoothing  the  image  just  along  radial  directions 
would  be  more  appropriate,  to  avoid  smoothing  along  wedges.  In  our  experiments  this  possible 
effect  was  not  noticed,  while  the  substantial  speed/efficiency  of  filtering  on  x  and  y  (and  not  along 
r)  encouraged  us  to  use  an  overall  Gaussian  blur. 

C.  Results  on  Real  Images 

The  threshold  for  deciding  N  of  the  junction,  rc,  was  set  to  0.4,  and  the  threshold  for  deciding  Rq 
and  ,  rw,  was  set  to  2.1,  for  all  the  images  we  used5.  We  observed  that  using  the  number  of 
discretization  units  of  the  angle  as  16  has  worked  well:  8  was  found  to  be  too  coarse  and  32  did  not 
result  in  any  significant  higher  accuracy6. 

Figure  8  shows  the  function/template  that  is  fit  to  different  points  on  the  images.  Figure  9  shows  an 
image  where  the  Adjunctions  have  been  detected  and  Figure  10  shows  the  results  of  Kona  on  a  region 
of  the  image.  After  the  filtering  it  computes  the  parameters  for  about  sixteen  junctions  in  a  minute, 
on  a  Sun  Sparc  Station. 

We  demonstrate  the  stability  of  the  junction  detector  by  running  it  on  a  pair  of  stereo  images  as 
shown  in  Figure  11.  (The  results  of  this  is  intended  for  use  by  stereo  algorithms.)  The  analysis  of  the 
result  of  this  is  presented  in  Figure  12.  We  measure  the  stability  by  first  finding  a  correspondence 
between  location  ( x\ ,  yf)  in  the  left  image  and  ( xr,  yr )  in  the  right  stereo  image  and,  then,  comparing 
the  number  of  junctions  in  a  small  neighborhood  of  the  two  locations. 

V.  Conclusions 

To  conclude,  we  summarize  the  contributions  of  this  paper: 

•  A  modeling  of  the  junction  that  includes  the  gradient  term  in  the  model,  and,  supports  the 
removal  of  a  small  disc  at  the  center  of  the  junction. 

5Recall  that  rc  and  rw  are  thresholds  for  relative  errors,  hence  are  unitless. 

6  We  choose  a  multiple  of  4  because  of  the  four  quadrants  -  using  any  arbitrary  number  could  do  as  well,  in  principle. 
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•  Using  the  model  to  obtain  the  scale  and  location  of  a  junction. 

•  A  simple  dynamic  programming  solution  to  the  optimization  problem. 

This  paper  demonstrates  the  successful  use  of  piecewise  constant  functions,  using  a  minimum  de¬ 
scription  length  principle,  to  detect  features  like  junctions  and  corners.  Similar  piecewise  constant 
functions  may  be  used  for  features  like  bar  detectors,  blobs,  end-points  etc. 

Kona  is  an  implementation  on  this  model.  This  detector  is  being  currently  tested  in  a  multitude  of 
tasks  such  as  object  recognition  (based  on  geometric  hashing),  detection  of  illusory  contours,  and 
depth  estimation  from  stereo  pairs  of  images.  It  is  also  being  tested  in  a  novel  approach  to  image 
compression  using  image  diffusion  and  multi-junctions. 
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Appendix 


I.  Computing  E 

Exploiting  the  radial  invariance  of  the  templates  T(0)  and  e(6)  we  can  write  equation  (7)  as  follows: 

POO  P  27T 


nZ7T 

[/2(r,  6)  -2  T(9)  I(r,6)  +  T2(6)]g(r)  rdrdO 

+xT  C[{-dnr,e]'i  °“'a ldI{r,e} 


de  f~um^FJ-  +  #mg%r)rdrde 


de 


=  C  CT\r,e)  + \(^p-)2]g{r)r  dr  de 


dl(r,0)  2. 
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N 

+  lim^o 
p= i 


2TpI(r,8)  +  \(re2(8)  -  2re(8) ^+2  —)\g(r)rdrdd, 


where  we  have  assumed  g*(r)  =  r2g(r)  We  can  now  exploit  the  radial  invariance  of  the  templates  to 
perform  the  integrals  over  dr.  Define  the  following: 


T 

C 

c * 

m 


fo  fQ  [l2(r,0)  +  \  (^A)  ]a(r)rdrd8, 


g(r)r  dr  =  Ri  -  R0 , 

i  foe  f? 3  _  d3 

g*(r)rdr  =  J  r2  dr  =  1  ^ . 

) 

/(r,  0)g(r)r  d  r. 


(13) 


Note  that 


rdl(r,e )  ,  ,  J  d  j'oo  QUO) 

Jo  ~r9{r)rdr=deJo  nr,0)9(r)rdr  =  —. 

The  radial  integrals  can  be  performed  numerically.  Thus, 


N 


E  =  T  +  Y  lira^o  f°P+'  f-2  TpI(0)  +  T2C  +  A(-2  re(8 )^~  +  C*e2(8))\d8. 


P  =  1 


We  can  now  perform  the  integrals  over  dO  to  obtain  the  following: 


N 

E  =  r+^2[(ep+1-ep)(-2TpIep,ep+l+T2C)  +  X(-2rep(deI)ep  +  C*e2p)l 

P  =  1 

where  we  have  used  that  e(0)  =  ep$(@  ~  0p) ,  and, 


(del)ep 

epF(0p) 

^  Op, i 


dl(°) ,  _  <9/0°°  l{r,0)g{r)rdr 


dO  1  p 

Op+i-e 


lim 

e->0 


{@p+i  Op)  J  Op 

3  P+l 

A8  Y  hi  AS). 


86 

e(8)F(8 )  d6 , 

Vn  , 

1(8)  de, 


3—3p 

where  the  integral  J^+1  1(0)  dO  is  adapted  to  the  image  lattice,  by  computing  1(0 )  at  every  interval 
AO  and  with  (0p+1  -  0p)  =  (jp+1  -  jp)A0. 
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II.  Numerical  Approximation  of  1(6) 

Using  the  definition  of  g(r)  =  1/r  and  equation  (11),  we  have 

I(°)  =  Jo  I(r,6)g(r)r dr 

=  i?/M)G)rdr 

To  carry  out  the  numerical  approximation  we  define  a  grid,  as  shown  in  Figure  13,  such  that  each 
cell  in  the  grid  corresponds  to  a  pixel  on  the  image.  A  ray  centered  at  (x0,  yo)  of  length  R\  is  broken 
by  the  grid  into  k  segments  !i,  . . ,  Ik  where  each  segment  lies  in  a  rectangular  cell  (pixel)  with 

intensity  ij,  j  =  1,2,...,  k.  Note  that  R\  =  /i  +  . .  .  +  4-  Thus  we  have  the  following  approximation: 

Further,  if  Rq  >  0,  the  segments  are  /m,  /m+i,  •  •  -  -Ik-,  f°r  some  1  <  m  <  k,  and,  R\  =  lm  +  /m+i  + 
•  •  •  +  %  and  I{6)  ss  Ei=m+i 

□ 
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Intensity 


*1  ^2 


Fig.  1.  Piecewise  constant  features.  A  junction  detector  is  shown  on  the  top  row:  3-junction  template  with 
edges  at  angles  9 i,  $2  and  63  with  intensities  I2  and  T3.  A  bar  detector  is  shown  in  the  bottom  row: 
it  is  characterized  by  an  angle  a  which  gives  the  inclination  of  the  bar  and  the  intensities  T\  and  T2  at 
“heights”  Y\  and  ¥2  from  the  center  of  the  location. 
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Marked  Images. 


er  versus  r  plots  for  the  three  locations. 


The  computed  windows. 


o 

Location  1  (46,50) 

$ 

Location  2  (68,81) 

O 

Location  3  (71,27) 

Rq  =  1,  R\  =  4  pixels. 

Rq  =  2,  Ri  =8  pixels. 

Rq  =  1,  R\  =  4  pixels. 

Fig.  2.  Estimating  Rq  and  Ri,  the  size  of  the  hole  and  the  size  of  the  window  respectively.  The  relative 
values  of  7 Z,  er  is  plotted  against  the  radius  r.  The  threshold  value  of  rw  =  2.1  was  used  to  obtain  the 
values  of  Rq  and  R\  as  shown.  Notice  that  when  we  pick  the  lower  bound  we  take  1  unit  less  than 

the  picked  interval  (for  instance,  the  range  for  the  first  graph  is  r  =  2  .  .  .4  and  Rq  =  1,  Ri  =  4).  We  do 
this  since  we  sample  the  values  of  77  only  at  integral  values  and  7lr  —  1ZS ,  r  —  1  <  s  <  r  where  s  is  a 
rational. 
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(a)  Position  on  the  image. 

(c)  Adaptive  window  at  (a) . 
(Ro  =  1,  Hi  =  3  pixels) 


(e)  The  best-fit  to  (c). 


(d)  Adaptive  window  at  (b). 
(Ro  =  2,  Ri  =  6  pixels) 


(f)  The  best-fit  to  (d). 


Fig.  3.  An  example  to  show  the  dynamic  computation  of  Rq  and  R\  at  different  locations  on  the  image.  We 
have  used  tw  =2.1. 


0 

The  junction  at  location 

Fig.  4.  Estimating  the  location  of  a  junction:  The  use  of  lZr ,  to  locate  the  center  of  the  junction.  L x>y 
indicates  the  x  and  y  coordinates  on  the  image.  For  convenience,  the  input  data  is  normalized,  that  is, 
the  intensities  are  scaled  to  lie  in  [0,  1].  The  unit  of  each  value  is  (intensity)2.  Note  that  the  location 
L/52,26  has  the  minimum  value  in  the  neighborhood.  [rw  =  2.1) 
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Fig.  5.  Estimating  the  number  of  wedges:  rc  is  used  to  threshold  the  relative  errors  to  declare  a  iV-j  unction. 
The  asterisk  denotes  the  computed  optimal  N  using  rc  =  0.4.  The  unit  of  E  is  (intensity) 2 ^length. 
The  thresholding  process  works  in  the  following  manner:  we  examine  the  successive  relative  errors,  in 
increasing  n,  until  we  reach  a  n  where  the  relative  error  exceeds  the  threshold,  rc.  The  previous  value  of 
n  is  declared  as  the  optimal  value. 
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Noise  in  image 


Fig.  6.  Test  of  stability:  The  range  of  intensities  in  the  image  is  0  —  255.  A  noise  with  a  Gaussian 
distribution,  with  standard  deviation  a,  is  introduced  randomly  at  each  pixel.  The  first  image  in  the 
top  left  has  a  =  0  (that  is,  no  noise)  and  the  angles  of  the  3-corner  at  the  center,  {x,y),  are  aij  =  90°, 
a2f  =  180°  and  =  315°  with  intensities  iij  =  120,  i2j  =  200  and  i3j  =  40.  a  is  varied  from 
0  to  69  to  obtain  24  images,  three  of  which  are  shown  in  the  top  row  for  illustration.  The  3-corners 
around  the  center,  (xa,ya),  of  each  image  is  extracted  and  the  following  errors  are  computed:  error  in 
location,  errLa  =  \xa  —  x\  +  \ya  —  y |,  error  in  intensities,  ( errla )2  =  ~  / ) 2 / 3 ,  error  in  angles, 

(errAfj)2  =  Y?j=i(ajo  —  cijj)2 f 3.  Finally,  errLa  vs  <r,  errAa  vs  a,  and,  errla  vs  a  are  plotted.  See  the 
text  for  further  explanations. 
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Fig.  7.  Test  of  stability  over  smoothed  noisy  images:  We  use  the  same  set  of  (noisy)  images  as  used  in 
Figure  6,  except  that  these  images  are  now  smoothed.  We  increase  the  smoothing  factor  as  the  image 
gets  noisier.  The  first  eight  images  (<r  =  0,  3,  . .  .,21)  are  smoothed  using  a  us  =  2.0,  the  next  eight 
(<r  —  24,  27, .  . . ,  45)  uses  crs  =  3.0  and  the  last  set  (<r  =  48,  51, ... ,  69)  uses  crs  =  4.0.  We  plot  the  errors 
as  in  Figure  6.  As  expected  we  are  able  to  get  rid  of  the  high  errors  for  the  very  sharp  images  (low  values 
of  <t). 
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(a)  Position  on  the  image. 

d 

(c)  Adaptive  window  for  (a). 
(Ro  =  1,  R\  =  4  pixels) 

J 

(e)  Projection  of  (c)  to  ID  signal. 

J 

(g)  Best-fit  to  (e). 


(b)  Position  on  the  image. 


(d)  Adaptive  window  for  (b) . 
(Ro  =  2,  R]  =  6  pixels) 


(f)  Projection  of  (d)  to  ID  signal. 


(h)  Best-fit  to  (f). 


Fig.  8.  Examples  of  a  sharp  corner  (column  1)  and  a  T-junction  on  an  image  (column  2).  We  have  used 
rw  =  2.1  to  get  Rq  and  R\  values,  and  rc  =  0.4  to  get  the  iV-j unctions. 
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(a)  (b)  (c)  (d) 


Fig.  9.  The  image  of  a  transparent  object  occluding  a  dark  object  in  the  background.  We  show  two  X- 
j unctions  detected  and  reconstructed  by  the  algorithm,  (a)  The  marked  images,  (b)  The  window  around 
the  marked  points,  (c)  The  corresponding  ID  signals,  (d)  The  4-junctions  at  the  points. 
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(a)  Marked  input  image. 


(b)  2-junctions.  (c)  3-junctions. 


Fig.  10.  Example  of  an  image  showing  the  results  of  the  junction  detector  Kona.  The  2-junctions  have  been 
filtered  which  has  removed  the  large  angle  corners  (  >  160  degrees).  Also,  the  low-contrast  3-junctions 
have  been  filtered  out.  The  values  of  the  various  control  parameters  are:  rc  =  0.4  (to  obtain  N  the 
number  of  wedges),  rw  =  2.1  (to  obtain  the  size  of  hole,  Rq,  and  the  radius,  Ri,  of  the  template),  number 
of  radial  lines  =  32. 
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(a)  Left  Image.  (b)  Right  Image. 


Fig.  11.  Kona  is  run  on  a  25  x  25  size  window  (at  location  (190,  350))  on  a  pair  of  stereo  images.  It  detects 
junctions,  with  significant  contrast,  on  similar  locations  on  both  the  images.  The  analysis  of  this  is 
presented  in  a  simplified  form  in  Figure  12. 
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Matrix  corresponding  to  the  left  stereo  image. 
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Matrix  corresponding  to  the  right  stereo  image. 

Fig.  12.  Kona  is  run  on  a  25  x  25  size  window  as  shown  in  Figure  11.  The  19  x  19  matrix  gives  a  count 
of  the  number  of  junctions  Kona  detected  in  a  7  x  7  window  centered  at  that  location.  (Thus  the  (1,  1) 
entry  of  the  matrix  corresponds  to  location  (193,353)  of  the  image  and  gives  the  number  of  junctions 
detected.)  By  inspection,  we  have  indicated  the  alignment  points  of  the  pair  of  images  by  presenting  the 
numbers  in  bold.  It  can  be  seen  that  the  number  of  junctions  (after  aligning  the  two  images),  match  up 
reasonably  well  (i.e. ,  gives  good  correlation  values). 
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Fig.  13.  The  pixels  shown  as  rectangular  cells.  The  ray  along  9  shown  by  the  solid  line  starting  at  (a?o,  yo)  is 
split  into  segments  |q  ,  l\, ...  ,ii.  Each  of  these  segment  lies  in  a  distinct  pixel  having  intensity  i$,i  i,  . .  . ,  U 
respectively. 


